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[57] ABSTRACT 

A multistage estimator is provided for the parameters of 
a received carrier signal possibly phase-modulated by 
unknown data and experiencing very high Doppler, 
Doppler rate, etc., as may arise, for example, in the case 
of Global Positioning Systems (GPS) where the signal 
parameters are directly related to the position, velocity 
and jerk of the GPS ground-based receiver. In a two- 
stage embodiment of the more general multistage 
scheme, the first stage, selected to be a modified least 
squares algorithm referred to as differential least 
squares (DLS), operates as a coarse estimator resulting 
in higher rms estimation errors but with a relatively 
small probability of the frequency estimation error ex- 
ceeding one-half of the sampling frequency, provides 
relatively coarse estimates of the frequency and its de- 
rivatives. The second stage of the estimator, an ex- 
tended Kalman filter (EKF), operates on the error sig- 
nal available from the first stage refining the overall 
estimates of the phase along with a more refined esti- 
mate of frequency as well and in the process also re- 
duces the number of cycle slips. 

5 Claims, 14 Drawing Sheets 
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MULTISTAGE ESTIMATION OF RECEIVED 
CARRIER SIGNAL PARAMETERS UNDER VERY 
HIGH DYNAMIC CONDITIONS OF THE 
RECEIVER 

ORIGIN OF THE INVENTION 

The invention described herein was made in the per- 
formance of work under a NASA contract, and is sub- 
ject to the provisions of Public Law 96-517 (35 USC 
202) in which the Contractor has elected not to retain 
title. 

1. Technical Field 

The invention relates to a multistage scheme for esti- 
mating the parameters of a received carrier signal, such 
as a carrier signal phase modulated by unknown data 
experiencing very high Doppler, Doppler rate, etc., 
that arises, for example, in the case of a Global Position- 
ing System (GPS) where the signal parameters are di- 
rected related to the position, velocity, acceleration and 
jerk of the GPS receiver relative to three or four satel- 
lites in the system. 

2. Background Art 

The problem of estimating the parameters of a re- 
ceived quasi-sinusoidal signal in the presence of noise 
occurs in diverse scientific and engineering disciplines. 
The signal parameters of interest are usually the phase, 
frequency and frequency derivatives which are varying 
with time. The estimation problem becomes considera- 
bly more difficult if the received carrier is modulated by 
unknown data while simultaneously experiencing very 
high Doppler and Doppler rate. As just noted, situa- 
tions occur in the cases of Global Positioning System 
(GPS) receivers, but may also occur in NASA deep 
space communication links under high spacecraft dy- 
namics. 

In a paper by W. J. Hurd, J. I. Statman and V. A. 
Vilnrotter tilted "High Dynamic GPS Receiver Using 
Maximum Likelihood Estimation and Frequency 
Tracking,” IEEE Trans., Vol. AES-23, No. 4, pp. 
425-437, July 1987, an estimator scheme is proposed 
and analyzed for the GPS system based on the maxi- 
mum likelihood estimation (MLE) of code delay and 
Doppler frequency over a single data bit period. This 
scheme estimates Doppler frequency (assumed constant 
over successive intervals of bit periods) and then deter- 
mines frequency and frequency rate by a Kalman filter 
tracking Doppler. The scheme does not involve carrier 
phase estimation. For the dynamic trajectories simu- 
lated in the paper by W. J. Hurd, et al., the approximate 
MLE performance exhibited a threshold of about 30 
dB-Hz in terms of the received carrier power-to-noise 
power spectral density ratio (CNR), below which rapid 
performance deterioration occurred. 

For GPS applications, an alternative scheme has been 
proposed by C. E. Hoefener and L. Wells in a paper 
titled “Utilizing GPS for ultra-High Dynamic Vehicle 
Tracking in Space,” Proceedings of the International 
Telemetering Conference, Las Vegas, pp. 771-773, 
October 1986, wherein a parallel (nondynamic) link is 
established between the GPS satellites and a control 
ground receiver for the purpose of communicating the 
data to the ground receiver. The ground receiver simul- 
taneously receives the frequency translated version of 
the GPS receiver signal and removes the data modula- 
tion from this dynamic signal. Such an effectively de- 
modulated signal is then processed by the estimation 
algorithm to obtain the desired signal parameter esti- 
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mates. There are several estimation schemes in the liter- 
ature for this problem. See, for example, R. Kumar, 
“Fast Frequency Acquisition via Adaptive Least 
Squares Algorithm,” Proceedings of the International 
5 Telemetering Conference, Las Vegas, pp. 91-101, Oc- 
tober 1986; R. Kumar, “Fast Frequency Acquisition via 
Adaptive Least Squares Algorithm,” IEE Proceedings, 
Vol. 136, Pt. F, No, 4, pp. 155-160, August 1989; R. 
Kumar, “Differentia] Sampling for Fast Frequency 
^ Acquisition via Adaptive Extended Least Squares Al- 
gorithm,” Proceedings of the International Telemeter- 
ing Conference, San Diego, pp. 191-201, October 1987. 

More recently in R. Kumar, “Efficient Detection and 
Signal Parameter Estimation with Applications to High 
Dynamic GPS Receiver,” JPL Publication 88-42, Na- 
tional Aeronautics and Space Administration, Jet Pro- 
pulsion Laboratory, California Institute of Technology, 
Pasadena, Calif., a scheme for simultaneous detection 
2 q and estimation has been proposed. This scheme is based 
upon first estimating the received signal’s local (data 
dependent) parameters over two consecutive bit peri- 
ods, followed by the detection of a possible jump in 
these parameters. The presence of the detected jump 
25 signifies a data transition which is then removed from 
the received signal. This effectively demodulated signal 
is then processed to provide the estimates of the global 
(data independent) parameters of the signal related to 
the position, velocity, etc., of the receiver. From the 
30 simultaneous, it is seen that the scheme offers very sig- 
nificant improvement in terms of the required CNR 
over the approximate MLE algorithm of W. J. Hurd, et 
al., (1987) cited above. A key feature of this scheme, 
which has a computational complexity of about three 
35 times that of a single extended Kalman filter, is that to a 
certain extent the data detection is independent of the 
acquisition of the phase or frequency of the received 
carrier signal in contrast to the conventional decision- 
directed phase-locked loop receivers in which the data 
40 detector is an integral part of the loop and depends upon 
the acquisition of the carrier phase and/or frequency. 
Thus, under low CNR and/or high dynamic conditions, 
the loop may not acquire lock or frequently lose it dur- 
ing tracking. 

45 An object of this invention is to provide an efficient 
method for estimating the parameters of a received 
carrier signal without undue computational complexity 
utilizing a multistage scheme. 

50 STATEMENT OF THE INVENTION 

In accordance with the present invention, the first 
stage of a multistage estimator operates as a coarse 
estimator resulting in higher rms estimation errors but 
55 with a relatively small probability of the frequency 
estimation error exceeding one half of the sampling 
frequency (an event termed cycle slip). The second 
stage of the multistage estimator operates on the error 
signal available from the first stage, refining the overall 
60 estimates, and in the process reducing the number of 
cycle slips. 

The first stage algorithm is preferably selected to be 
a modified least squares algorithm referred to as differ- 
ential least squares (DLS) algorithm. This estimation 
65 stage provides relatively coarse estimates of the re- 
ceived signal frequency and its derivatives. The second 
stage algorithm is preferably a third-order extended 
Kalman filter (EKF) which yields a more refined esti- 
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mate of frequency as well as an estimate of the signal 
phase. 

A major advantage of the proposed multistage esti- 
mator is a reduction in the threshold on received carrier 
power-to-noise power spectral density ration (CNR) as 5 
compared to the threshold achievable by either of the 
two cascaded estimator stages alone. In fact, it appears 
from simulations that for the case of an unmodulated 
carrier, the proposed scheme achieves the same thresh- 
old as for an almost exact and computationally intensive 10 
implementation of the maximum likelihood estimator 
(MLE). For the case of a data modulated carrier, the 
proposed scheme provides an improvement of about 6 
dB in terms of CNR compared to an earlier approximate 
MLE scheme reported by W. J. Hurd, et al., (1987) 15 
cited above. The over-all complexity of the algorithm is 
about two times the complexity of a third-order Kalman 
filter or a single fourth-order extended Kalman filter. 

The novel features that are considered characteristic 
of this invention are set forth with particularity in the 20 
appended claims. The invention will best be understood 
from the following description when read in connection 
with the accompanying drawings. 

BRIEF DESCRIPTION OF THE DRAWINGS 

FIG. 1 illustrates a block diagram of a prior-art sys- 
tem for signal parameters estimation by differential least 
squares (DLS) algorithm. 

FIG. 2(a) illustrates a block diagram for signal param- 
eters estimation by least square (LS) algorithm with 
model noise colored. 

FIG. 2(b) illustrates a conceptual block diagram for 
signal parameters estimation by least square (LS) algo- 
rithm with prewhitened noise. 

FIG. 2(c) illustrates a block diagram for signal param- 35 
eters estimation by least square in an approximate realiz- 
able equivalent of FIG. 2(b). 

FIG. 3(a) illustrates a generalized single-stage estima- 
tor used in the first stage of the present invention. 

FIG. 3(b) illustrates a functional block model for a 40 
more compact representation of the generalized single- 
stage estimator shown in FIG. 3(a). 

FIG. 3(c) illustrates a functional block model for 
subsequent stages of a multistage estimator. 

FIG. 4 illustrates a multistage estimator for the pro- 45 
cess 6j(t). 

FIG. 5 illustrates in graphs a, b and c, a high dynamic 
trajectory used in simulated analysis of tracking perfor- 
mance of a multistage estimator shown in FIG. 4 
wherein n=2. 

FIG. 6 is a graph of the probability of losing fre- 
quency lock vs CNR for a DLS algorithm in the ab- 
sence of data modulation. 

FIG. 7 is a graph of RMS frequency estimation error 
vs CNR for a DLS algorithm in the absence of data 55 
modulation. 

FIG. 8 is a graph of the probability of losing fre- 
quency lock vs CNR for a DLS algorithm in the pres- 
ence of data modulation. 

FIG. 9 is a graph of RMS frequency estimation error 60 
vs. CNR for a DLS algorithm in the presence of data 
modulation. 

FIG. 10 is a graph of the probability of losing fre- 
quency lock vs CNR for a DLS-EKF algorithm in the 
absence of data modulation. 

FIG. 11 is a graph of RMS frequency estimation 
error vs CNR for a DLS-EKF algorithm in the absence 
of data modulation. 
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FIG. 12 is a graph of modulo 2 t t RMS phase error vs 
CNR for a DLS-EKF algorithm in the absence of data 
modulation. 

FIG. 13 is a graph of the probability of losing fre- 
quency lock vs CNR for a DLS-EKF algorithm in the 
presence of data modulation. 

FIG. 14 is a graph of RMS frequency estimation 
error vs CNR for a DLS-EKF algorithm in the pres- 
ence of data modulation. 

FIG. 15 is a graph of RMS pseudo range estimation 
error vs CNR for a DLS-EKF algorithm (with and 
without data modulation). 

DETAILED DESCRIPTION OF THE 
INVENTION 

1. Introduction 

An alternative scheme for the estimation of the signal 
parameters will now be described for both the case of 
unmodulated carrier signal and the case of a carrier 
signal phase modulated by unknown data. The pro- 
posed scheme is somewhat simpler than that of R. Ku- 
mar, described in the aforesaid JPL Publication 88-42, 
in that it is not essential to explicitly detect the data 
modulation for the second case. Basically, the new algo- 
25 rithm involves an appropriate modification of the DLS 
scheme of R. Kumar in the aforesaid paper titled “Dif- 
ferential Sampling for Fast Frequency Acquisition via 
Adaptive Extended Least Squares Algorithm” so as to 
apply the algorithm to the case of unknown data modu- 
30 lation. 

As discussed in that paper, if a DLS technique is 
applied with the Nyquist sampling of the received sig- 
nal, a loss in performance is expected compared to the 
optimum achievable performance. Consequently, over- 
sampling and cyclic sampling was proposed to avoid 
such a loss. In the present invention, sampling at Ny- 
quist rate is proposed with an alternative method for 
estimating parameters in order to keep the overall per- 
formance close to optimum. The new scheme proposed 
consists of a multistage procedure wherein the parame- 
ters of the signal are estimated in more than one stage. 
First, coarse parameters are estimated by an algorithm 
like DLS which has a low threshold on CNR but with 
possibly higher rms estimation errors. Then an error 
signal whose parameters are equal to the difference 
between the true parameters and the coarse estimates is 
processed by another algorithm to estimate these error 
signal parameters. Since the error signal involves much 
smaller dynamics, the second algorithm can have 
50 smaller bandwidth resulting in a smaller estimation 
error. In principle, this procedure of processing an error 
signal in another stage can be repeated any number of 
times with successive stages having progressively lower 
bandwidths. 

The example described below of a multistage estima- 
tor is confined to two stages of recursion and applies a 
modified least square algorithm for the first stage and a 
third-order extended Kalman filter algorithm for the 
second stage. It was expected that the overall algorithm 
would have both smaller threshold and smaller estima- 
tion errors compared to either algorithm operating by 
itself. Indeed, this is borne out by simulations presented 
hereinafter. Thus, for the case of no data modulation, 
the threshold on SNR is about 1.5 dB lower than the 
65 third-order EKF, and the estimation errors are only 
marginally higher than for the third-order EKF alone. 
The threshold achieved is in fact the same as achieved 
for a nearly exact implementation of the maximum like- 
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Iihood estimator (MLE). It is also noted that the thresh- 
old achieved by the present invention is about 3 dB 
lower than conventional cross product AFC (CPAFC) 
loops and phase-locked loops compared by V. A. Viln- 
rotter, S. Hinedi and R. Kumar, “A Comparison of 
Frequency Estimation Techniques for High Dynamic 
Trajectories,” JPL 88-19, Jet Propulsion Laboratory, 
California Institute of Technology, Pasadena, Calif., 
whereas the rms error is less than one half of that ob- 
tained by CPAFC. In the simulation described hereinaf- 
ter, the rms error is marginally higher than for a third- 
order EKF due to the nonoptimal sampling used in the 
DLS algorithm. 

For the case of data modulation, results of the simula- 
tions are compared with those reported by Hurd, et a!., 
(1987), where analysis and simulations are presented on 
the performance of Fast Fourier Transform (FFT) 
based MLE algorithm. In that report, the trajectories of 
the GPS signals have somewhat less severe dynamics 
compared to those considered in this presentation. In 
terms of threshold on CNR, the proposed scheme of 
this presentation exhibits a threshold of 24 db-Hz com- 
pared to about 30 db-Hz reported by Hurd, et al., thus 
providing an improvement of about 6 dB. In terms of 
the rms frequency estimation errors at a 30 db-Hz CNR, 
the scheme of Hurd, et al., provides a rms range rate 
error of about 6 m/s compared to an error of less than 
2m/s achieved in this simulation. There is also a very 
significant improvement in terms of the rms position 
estimation error. At about 30 dB-Hz, an rms error of 1 
meter is reported by Hurd, et al., compared to about 
0.25 meter obtained by the proposed algorithm. It may 
also be remarked that in the previous scheme of Hurd, 
et al., pseudo-random codes with rate 10.23 MHz are 
needed for the purpose of range measurements, thus 
requiring a zero-crossing channel bandwidth in excess 
of 20 MHz. The present invention on the other hand, 
extracts the range information from the carrier signal 
itself and thus needs a bandwidth equal to only a frac- 
tion of 1 MHz. 

2. Receiver Configuration 

Consider the problem of estimating the high dynamic 
phase process 0Xt) of the desired signal sXO observed in 
the presence of an additive narrow-band noise process 
vXt) as 


of a signal sz,(t) from a voltage-controlled oscillator 
(VCO) 14 and lowpass filters IS and 16 as shown in 
FIG. 1. The VCO frequency update rate of 1/T S indi- 
cated by a block 17 is controlled by a signal /J(N) pro- 
5 duced by a differential least square (DLS) algorithm 
indicated by a block 18 the function of which will be 
described below. Assuming that the input of the VCO is 
a signal 6l( t) which is an appropriate quadratic function 
of time, t, resulting in the following VCO output signal 
10 Si (t) 


rKO = + >7(0 

= /ISintacl + 0/(/) + ttD(i)) + vj(t). 


( 1 ) 


Sift) =2 Cos (o> c t+0L(t)) 

Sift) = 0£ o +M£ o +7Lot 2 + S^t 3 


(3) 


15 


for some constant vector \pLo=l&Lo o>Lo JLo 8l 0 ]', the 
sampled version of the in-phase and quadrature compo- 
nents of the demodulated signal are given by 

20 y(k)=A Sin (0(k)+wD(k))+v,<k) 

z(k)= ACOS(0(K) + irD(K)) •+■ v^k); k=l, 2 N (4) 


where 


25 0(k) = 0/k) - «L(k) =e 0 + ct>okT,-f 7o(kT,) 2 ‘ 

+ 5„(kT,) 3 

'l>o=4'/o->l>Lo=[0o °>o Vo So]' 

30 and i ji 0 is the parameter vector characterizing the error 
signal to be estimated, with T s denoting the sampling 
interval. In Equation (4), v,(k) and v ? (k) represent the 
sampled in-phase and quadrature components of the 
bandpass noise process vXk). The parameter vector is 
35 estimated by the differential least squares (DLS) algo- 
rithm of R. Kumar (1988) cited above, as described in 
the following section. 

3. Differential Least Squares (DLS) Algorithm 
Consider first the problem of estimating the unknown 
40 parameters a> 0 , 7o and d 0 from the measurement Equa- 
tion (4) for the case of no data modulation, i.e., when 
D(k)=0, and expanded Sin(0(t)) in a Taylor series 
around U--i =(k— 1)T 5 to obtain 
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Sin(0(t))=Sin (0(k- 1)+T0(k- 1) Cos (0(k- 1))+. . 


(5) 


where a> c is the received signal carrier frequency in the 
absence of any dynamics and D(t) is a binary digital 
waveform. In the case of a GPS receiver, the dynamic 
phase process 6X0 arises from the receiver dynamics 
and, over a sufficiently small estimation period, 


with a similar expansion for Cos (0(t)). For small (t-t*-i), 
the series in Equation (5) may be approximated by the 
50 first two terms and from Equation (4) the following 
differential signal model with T*=(k-($))T*: 

A 

ya(k) = y(k) - y (k - 1 ) = TXcflo + Zvor* + - 
38^* 2 )z(k-l) + f,<k) 
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0X0 = 0/ o + <u/ 0 t 1 + y / 0 t 2 4- 8/ot 3 


z^k)=z(k)-z(k-l)=Ts(£Oo+27or/r+38oT*. 

( 2 ) 2 )y(k-l)+?Xk) 


( 6 ) 


for some unknown parameter vector \fii 0 =[&io vio 77c 
S/ o], In a somewhat simpler version of the problem, the 
data modulation D(t) is either absent or is assumed 60 
known and thus can be eliminated from Equation (1). 
Both of these cases will be treated in some detail. 

In the first stage, the present invention estimates the 
parameters related to the frequency and its derivatives 
using the DLS algorithm. For this purpose, the re- 65 
ceived signal sXt) is quadrature demodulated in a con- 
ventional input section 10 comprised of mixers 11 and 
12, a tt/ 2 phase delay element 13 which shifts the phase 


where 

liW=v<k)— v/k-l)— Tj(<i> 0 + 27 0 r*+- 
38 0 Tt 2 )v 9 (k— 1) 

i(k) = v^(k) - v^(k - 1) -Tj(a) 0 + 27 0 rk+ 3 S 0 Tk 2 ' 

)v,(k — 1) (7) 

The measurement model of Equation (6) may be rewrit- 
ten in the following standard: 

Zd(k) = H'(k)/S + €(k); k = 1 . 2, . . . N (8) 



where 


)(A 2 + cr v 2 ) l [\TkTk 2 \'T7 X 


(14) 


P' = [a 0 2y 0 68„] 


5 


H(k) = 


TAk - 1) 
-Ts)Kk - 1) 


T s r k z(k - 1) 0.5T s T k Mk - 1) 
-T s Tky(k - 1) -0.5T s Tk 2 y(k - 1) 


(9) 


10 


Z' < Xk)=[y^k) Zd(k)]; 4'(k)=[|Xk) «#)] 
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The parameter vector P in Equation (8) is now esti- 
mated by an exponential data-weighted least squares 
algorithm in a recursive or nonrecursive from (Kalmal 
filter). In its nonrecursive form, the estimate of P ob- 

20 

tained on the basis of N measurements and denoted by 
P(N) is given by 


( 10 ) 


/3(A) 


= ^ # mj)lf(j)\ N -j 'j [j‘=\ H <f> z ^ xN ~ J } 


25 


where X is some appropriate weighting coefficient 
within 0<X< 1. An equivalent recursive form of Equa- 
tion (10) is the following algorithm: 


/3(N)=/3(k-l)+L(k)«(k) 

L(k) = P(k - l)H(k)[XI + H'(k)P(k — 1) H(k)]‘ 1 

P(k) = {P(k - 1) - P(k - 1 )H(k)[XI + H'(k)P(k -1)H(- 
(k)]-'H'(k)P(k+l)}/\ 


35 


OD 


€(k)=Zz(k)— H'(k)fl(k — 1); k=l, 2 N 


40 


Note that the matrix to be inverted in Equation (11) is 
only a (2x2) matrix. In an alternative but equivalent 
form, one may process the scalar measurements y^fk), 
Zd(k) sequentially instead of working with the vector 45 
measurement Zj(k). Moreover, the matrix P(k) of Equa- 
tion (1 1) with k=N is the same as the matrix inverse in 
Equation (10), i.e., 


f-'(k) 


■- 1 H(jW(j)k k -J 

j-1 


( 12 ) 
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Alternatively, the matrix P-'(k) may be written as 


55 


k A 

P-'(k) = \ k -j{^(/) + /Q^BjTs 2 ; Bj = 


r- T 2 r 2 

T J T 1 T J 


T 2 t 3 t 4 
T J T J T J 


(13) 


60 


and thus the matrix P(k)H(k) required in the update of 
/3(k), and equal to L(k) in Equation (1 1), may be approx- 
imated by 

P(k)H(k)=X(k)[z(k)-y(k)] 


X'(k) ■■ 




where 

E[v,2]=E[v ? 2]=(i)<r v 2; E[z2(j)+y 2 ©] = A 2 + <r v 2 . 

In Equation (14) the vector X(k) is data independent 
and thus could be precomputed for k=l, 2, . . . , N for 
computational simplification of Equation (11). Simi- 
larly, in the implementation of Equation (10), the first 
matrix may be replaced by the data independent matrix 


^ + o.2)T 2 . 

3(a). Modified Least Squares Algorithm 

It the noise £(k) in the signal^ model of Equation (8) 
were white, then the estimate /3(N) obtained from the 
algorithm (10) or (11) would approach 0 as N— ►*, if 
one ignores the approximation made in arriving at 
Equation (8) and X is selected equal to 1. However, as 
the noise £j(k) in Equation (8) is colored, there would be 
considerable bias in the parameter estimates under low 
to medium signal-to-noise ratios. To reduce such a bias 
or possibly eliminate it, the following simple modifica- 
tion is proposed. If the instantaneous frequency co(t£) at 
time n given by (w 0 +2y 0 Ti+3S 0 T/ : 2 ) appearing in 
Equation (7) is small compared to \/T s , then the noise 
vector £(k) is equal to v(k) — v(k— 1) where v (k)=[v,(k) 
v ? (k)]'. This situation is illustrated in FIG. 2(a). 

To eliminate the bias, the noise £(k) must be whitened 
by passing through the transfer function (1— z -1 )- 1 as 
shown in FIG. 2(b) where Pub denotes an unbiased 
estimate of p. The least squares algorithm is, in general, 
nonlinear and time-varying. However, if it is assumed 
that the algorithm in Equation (11) asymptotically ap- 
proaches a time-invariant system, then under such an 
assumption, one may interchange the least squares algo- 
rithm with the transfer function (1— z -1 ) -1 to arrive at 
the arrangement of FIG. 2(c). This, of course, corre- 
sponds to post-averaging the least squares estimates, 
which is the preferred technique for the present inven- 
tion. Such a simple procedure provides very significant 
improvement in the estimates P( k) when the signal-to- 
noise ratio is low. In the simulations of the next section, 
the infinite time average (1— z -1 ) -1 of FIG. 2(c) is 
replaced by an exponentially data-weighted averaging 
to take into account the time variation of the parameters 
to be estimated. 

3(b). DLS Algorithm in the Data Modulation Case 

In this case, the data samples D(k) in the signal model 
of Equation (4) take only possible values ± 1 and the 
received signal may equivalently be written as: 

y(k)=D(k)A Sin (fl(k))+vXk) 


z(k) = D(k)A Cos (0(k)) + V? (k); k = 1 , 2, N (15) 

Thus, as may be easily verified over any bit interval T*, 
where D(k) remains constant, the differential signal 
model of Equation (6) remains valid irrespective of the 
value of D(k). The model of Equation (6), however, is 
not applicable for those samples which lie on the bit 
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boundaries, i.e., when y(k) and y(k— 1) lie in different 
bit intervals. A simple modification of the algorithm to 
take care of the data modulation case is to simply dis- 
card such differential samples. If the number of samples 
M over any bit period is fairly large, this would incur a 
negligible loss in the effective signal-to-noise ratio com- 
pared to the case of no data modulation. In fact, such a 
loss is simply equal to 10 logio(l-l/M)dB which is 0.45 
dB for M = 10. This is corroborated by the simulations 
of the next section. 

• • 3(c) Estimation of Time-Varying Parameters 

In the signal model considered above, it is assumed 
that the input signal parameter vector i p/o is either a 
constant or a slowly varying function of time. In prac- 
tice, this may be the case only over relatively short 
intervals of time, but there may be large variations in 
vlf/oover a comparatively large observation period. To 
take into account such a variation and to ensure that the 
instantaneous difference frequency fi(t) = d/dt(0(t)) 
(the sampled version of 0(t) given in Equation (4)) re- 
mains within the low-pass filter-pass band of FIG. 1, the 
parameter vector i J//, 0 generating the instantaneous fre- 
quency of the VCO is updated at regular intervals of 
T=NTjsec for some integer N. The parameter vectors 
tpLo and >!> 0 would change their values at intervals of 
T sec, assuming that the value of N is selected to be 
sufficiently small so that the variation in i jii over any T 
sec interval is small. Denoting by &Lo(T +),a>Lo(J +), 
etc., the values of reference oscillator parameters just 
after the update at time T, are: 

0l<XT+)=0z.o(T-) 

oilJJ + =fli«(T— )&a(0/T)+2Ty^0/T)-t-3T- 
2 8„(0/T) 

7 Lo(J + ) = rio(T - ) + 3T5o(0/D 

8£o(T+) = A £o (T-) + 8 o (0/T) (16) 

In Equation (16), 0[_ o (T—), etc., represent 

the oscillator instantaneous phase, frequency, etc., just 
before the correction, and the remaining terms on the 
right hand side represent the correction made on the 
basis of the estimation algorithm. Thus, 

0L O (T - e Lo (P+ ) + 0 >lMP + )T + 7 L„( 0 + )T 2 +- 

8£o(0+)T 3 

n lo( r - ) = aiAO + ) + 2y lo(0 + )T + 38l„(0 + )T 2 
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z(N+j)=A Cos (0(N + J)) + v,(N+j); j = l N 

<?(N+j) = 0„(T+O 

l+UotT+hT^+yotT+XjT^+Soa+XjT,) 2 (18) 

^ The last three elements of the vector \p a (T+) will be 
zero if there is no change in the input signal parameters 
over the T sec interval and the estimate of <j» o (0+) is 
obtained with zero estimation error. Thus, the tech- 
nique is to set the a-priori estimate of the vector i|/ 0 (T +) 

10 equal to 0 and apply the DLS algorithm to estimate 
ip 0 (T+) on the basis of observations fy(N+j), z(N+j); 
j = l, . . . , N}. The measurement model is obtained by 
simply replacing the index k by k+N in y(k),z(k)- 
in Equations (6-9) but with n— (k-ji)T f as 
before (corresponding to a shift in time reference). 

In' the estimation of i|/ 0 (T +) via the recursive algo- 
rithm of Equation (1 1) with the index k=N + 1, . . . 2N, 
the “initial” covariance matrix P(N -f 1) is obtained as: 

20 P(N + 1) = XFP(N)F' + Q (19) 

where 


A 1 2 T 3 T 2 
F= 0 1 IT 
0 1 1 

30 and the matrix Q represents the uncertainly introduced 
due to the change in the input process parameters over 
the interval of T sec. Specifically, the last diagonal 
element of Q represents the variance of the change in 
the parameter 68 (equal to the second derivative of 
35 frequency and related to the jerk of the physical trajec- 
tory) over the interval T. The above procedure is then 
extended in a straightforward manner to the subsequent 
update intervals. The estimates of the input signal phase 
and frequency at time instances 1T+ are then simply 
40 given by 6i 0 ( IT +) and o>lo( IT +) respectively for 1=0, 
1, 2 

4. Multistage Estimation 

Most of the phase and frequency estimation schemes 
for the first stage of a multistage estimate can be repre- 
45 sented as in FIG. 3(a) in which elements common with 
FIG. 1 are identified by the same reference numerals. 
An update algorithm (represented by a functional block 
21) for generating the correction signal for the VCO 
could be any recursive or semi-recursive algorithm 


r(T-)=y Lo (04-)+3T8 £ ^0-(-) 

AUT- )=8 lo (0 + ) (17) 

Note that in Equation (16) there is no correction in the 
oscillator phase as the DLS algorithm does not provide 
the phase estimate. In Equation (16) w o (0+T) denotes 
the estimate of parameter a> 0 ( 0+) obtained on the basis 
of measurements up to time T. Since there is no step 
change in the oscillator phase, the sampled measure- 
ments y(N) and z(N) at the demodulator output are the 
same with or without a correction at the instance NT S . 
However, the subsequent measurements y(N+j) and 
z(N+j) are now expressed with respect to the new 
parameter vector &) 0 (T+) = [# 0 (T+) co 0 (T+) 8 0 (T+)]'- 
— ^/o(T+)— il»z.o(T+) as in Equation (18) below. 

y(N+j) = A Sin (0(N+J))+vKN+j) 


50 including an EKF or DLS algorithm, and the VCO 
update interval may be some integer multiple of the 
sampling period T s . FIG. 3(b) illustrates a functional 
block model for the generalized first-stage estimator 
illustrated in FIG. 3(a). This functional block model 
55 will be used hereinafter as a symbol to illustrate the first 
stage in a multistage estimator shown in FIG. 4 in ac- 
cordance with the objective of this invention. Note that 
0ft), v,(t) and v^ft) of FIG. 3(a) are respectively equal to 
0i(t), v,'(t), v ? i(t) of FIG. 3(8). FIG. 3(c) illustrates a 
60 functional block model for subsequent stages of the 
multistage estimator of FIG. 4 . 

An important observation made from FIG. 3(8) is 
that along with the estimate 0X l ) of the input phase 
process 0X0, there is also available a pair of signals, 
65 y'(t)= A Sin 041(t) + v, 1 (t), z'(t)=A Cos 0(t)-t-v 9 '(t) 
dependent upon the estimation error 0(0 = 0X0 — #£(0- 
These error signals have exactly the same form as the 
signals at the input to the estimator. Moreover, the 
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additive noise associated with these signals has statistics 
identical with the statistics of the noise at the input to 
the estimator. Therefore, this leads to the interesting 
possibility of the present invention, namely estimating 
the error signals &(t) in a way similar to the estimation 5 
of 0/(t). In fact, the procedure can be repeated any 
number of times as shown in FIG. 4 , although only a 
two-stage estimator is described below as an example of 
the present invention. It should also be noted that only 
the first estimation stage shown in FIG. 3(a) requires a 10 
VCO for down conversion because the actual input to 
this stage is at rf frequency ia c , whereas subsequent 
estimator stages generate the error signals by simple 
baseband computations. For example, in the discrete- 
time version of the estimation procedure, the signal at 15 


the output of estimation stage m may simply be com- 
puted as 

y m( k) = y m-l( k ) Cos Sin 

<«"-'(k-l) 20 

z m (k)=y m -'(k) Sin (0 m -'(k-l)+z m “ 1 (k) Cos 
(9 m — '(k— 1); m = 2, 3 n (20) 

The refined estimate of the phase process 6/(t) in FIG. 

4 at sampling instance k is then simply given by 25 

exk)=« i (k)+e 1 (k)+ • • • +nk) < 2! ) 


. The advantage of such a recursive estimation procedure 
is that the overall threshold in terms of the required 30 
CNR for the multistage estimator can be made much 
smaller than a single stage estimator, especially in situa- 
tions involving high dynamics. 

In the prior-art single-stage estimator shown in FIG. 

1, the process parameters may be assumed to remain 35 
constant only over short intervals of time due to the 
high dynamics involved in a receiver carried, for exam- 
ple, in a high performance aircraft or spacecraft. Thus, 
the estimator is forced to use a relatively large noise 
bandwidth (shorter averaging period), resulting in large 40 
errors in the phase and/or frequency estimates. If the 
estimation errors are outside the region over which the 
error model (linear) assumed for the estimator remains 
valid, the estimator is said to be working below thresh- 
old or in the out-of-lock condition. In this condition, the 45 
estimation errors can be several orders of magnitude 
higher compared to the operation above threshold. 

In the multistage estimator of the present invention, 
this difficulty can be circumvented by successive reduc- 
tion of the dynamics (the estimation errors due to dy- 50 
namics) at the output of consecutive estimator stage and 
by averaging the signal over progressively longer inter- 
vals (and thus progressively reducing the effect of 
noise) over which the process parameters remain nearly 
constant. In this multistage estimation structure, none of 55 
the individual stages (except the last one) need necessar- 
ily operate above its threshold. For the convergence of 
the overall estimator, one only requires that the esti- 
mates are made in the right direction (estimation errors 
do not exceed the parameters to be estimated in some 60 
average sense). 

5. Recursive DLS-EKF Algorithm 

The following section considers a simple special case 
of n=2 for the more general multistage estimator of 
FIG. 4, wherein the first stage is as described above 65 
with reference to FIGS. 3a and 3b and the DLS algo- 
rithm is as described in the previous section above. The 
second stage is an extended Kalman filter (EKF) as 


12 

described above with reference to FIG. 3c. As the dy- 
namic variation of the error signal $(k)=0 ! (k) of the 
first stage is much smaller compared to the original 
signal 9( k) (the frequency variation over any update 
interval is much smaller), the effective averaging time 
period for the second-stage estimator implemented with 
a Kalman filter can be selected to be higher than for the 
first stage implemented with a DLS algorithm. This is 
achieved by selecting a smaller value of the “dynamic 
noise” covariance matrix Q for the EKF and/or a 
higher value for the exponential data weighting coeffi- 
cient X of Equation (9). 

5a. Extended Kalman Filter (EKF) Algorithm 
Referring to FIG. 4, consider the problem of estimat- 
ing the unknown error signal parameters a> 0 /, y 0 /and 8 0 l 
in the 1th VCO update period for any integer 1 g 1 on the 
basis of the set of measurements {y'(k),z'(k)} of Equa- 
tion (22) below by an EKF 

y’(k)=A Sin (e«(k))+v/>(k) 
z'(k)=A Cos (0'(k))+v ? ’(k) 

9'Qt)=0ol+ <*> 0 /(jT2) 2 +8oXjTj) 3 

k=N(l — l)+j; j = l, 2 N; 1 = 1, 2, . . . (22) 

Note that as for the first-stage DLS algorithm the pa- 
rameter vector i l) 0 i=[9 0 i <*oi Jot So/]' may be different 
over different VCO update intervals. For computa- 
tional simplicity, a third-order EKF is used in the sec- 
ond stage and the contribution of the last term in the 
expression for 0*(k) is ignored, which is appropriate for 
the GPS trajectories considered here. Denoting the 
state and parameter vectors at time k=N(l— l)+j by 
ii/Xj) and n/ resnerfi'"* 1 ” : with ij/(j)=[l jT* 

O.5(jT s ) z ]',i7; = [0 O (a)o;2Xo;]', the extended Kalman filter 
equations for the update of T); in the second stage, the 
estimate of t )/, are given by 

’)/ = ’!/(/- 1) + ( 23 > 
m = ikj - do.(/kx + vow - d«-‘ 

2/0) = (2/0 -1) - 2/0 - W(/)[X + 

fO)2io - nwir 1 'i> , 0)2/0 - ni/x 

vj(/) = y\k)Cos(<j>l(J)) _ z>(k)Sin(<£/(/)) 

ifli = 4V)viO' - 1 );* = (/— 1 W + J: 

j = 1,2 N; l = 1,2, . . . 

In the Equation (23) above, the initial estimate 77/(0) is 
simply taken to be equal to t)m(N). This is an appropri- 
ate choice for the initial estimate in view of the fact that 
if the first stage of the DLS estimation algorithm is 
covergent then tjj would possess some continuous drift 
term, i.e., w 0 i will have a component linear in time if 
. j The “initial error covariance” matrix 2/0) is 
simply set equal to some diagonal matrix representing 
the uncertainty in the difference parameter 
5b. Estimation in the Presence of Data Modulation 
In this case one could apply the ore sophisticated 
version of R. Kumar (1988) cited above, wherein an 
explicit detection of possible data transitions is followed 
by the demodulation of data, thus effectively reducing 
the problem to the case of no data modulation consid- 
ered above. However, here such a detection is by- 
passed, and instead a simple modification in the estima- 
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tion algorithm is proposed that takes into account the 
data modulation. If the VCO update interval T is se- 
lected equal to bit period T*, then the data modulation 
represents an additional phase uncertainty at the bound- 
aries of the update intervals. This is taken into account 5 
by adding an appropriate value, say (tr/2) 2 , to the first 
diagonal element of X/(0) and modifying the initial esti- 
mator f 7/(0) by 7 t/ 2, i.e., T}i(0)=7ji i(N)+7r/2, for those 
values of 1 that correspond to bit boundaries. Such an 
algorithm is expected to result in somewhat higher io 
•estimation errors compared to the more sophisticated 
scheme of R. Kumar (1988) cited above, but is much 
simpler in terms of implementation. In the case of a 
two-stage estimator (stages 1 and 2 of FIG. 4) , the 
estimates of the input signal phase and frequency at time 15 
instance IT are given by 

l(lT)=»i((I— l)T +)+<Ji’(N)iji(N) 

Q(lT)=ft£((l-l)T+)+S 0 ;+2w(NT J ) (24) ^ 

where oi 0 r and y 0 i represent the second and third ele- 
ment respectively of i)i(N)- 

6. Simulations 

Simulation results obtained when the algorithm is ^ 
applied to the tracking of phase and frequency for high 
dynamic GPS receivers will now be presented. For the 
purpose of simulation, assume that the pseudo-random 
code has been removed from the received signal, and 
symbol timing has been acquired. The simulation con- 
siders both the case of the data modulation removed via 3 
an auxiliary link and the case of an unknown modula- 
tion present. For the simulations a sampling rate of 500 
samples/second is assumed with a high dynamic trajec- 
tory considered previously in V. A. Vilnrotter, et al, 

“A Comparison of Frequency Estimation Techniques 3 
for High Dynamic Trajectories,” JPL 88-21, Jet Pro- 
pulsion Laboratory, California Institute of Technology, 
Pasadena, Calif., September 1988, and reproduced in 
FIG. 5 where the acceleration and the jerk (the deriva- 
tive of acceleration) are measured in units of g (the 
gravitational constant equal to 9.8 m/s). In the case 
when data modulation is present, a BPSK modulation at 
a rate of 50 bits/second is assumed. 

The parameters of most interest in this application are ^ 
the instantaneous phase and frequency of the input sig- 
nal r/t), which corresponds to the high dynamic GPS 
trajectory of FIG. 5. Since the present invention is 
mainly interested in the tracking performance of the 
proposed algorithm, it is assumed that the initial trajec- 
tory parameters at zero time are known as in V. Vilnrot- 
ter, et al., Sept. 1988, supra. The received signal carrier 
frequency fe=<D c /2jr in the signal model of Equation (1) 
is taken to be equal to 1.575 GHz. The GPS receiver 
instantaneous pseudo-rand R in meters and velocity Vd 
in m/s are related to the instantaneous phase drft) of " 
Equation (1) and its derivative 6j(t) as 


0 c_ 

2tt f e 


d ~ Jd fc ~ 2* fc 

where fe denotes the_Doppler_in Hz and c is the speed of 
light. Denoting by R(l) and fe(l) the estimates for R(l) 65 
and fe(l) respectively, which denote the range and Dop- 
pler of the input trajectory at the end of the 1th update 
interval, then the performance measures of the estima- 


tion algorithm are given by the following sample rms 
values of the estimation errors 


Xj im - R(t )] 2 


firms =77 \J / 2 [ Vd!) -fm 1 


where L=4000/N is the number of update intervals for 
the entire trajectory. These measures are obtained as a 
function of P/N 0 , where P denotes the received carrier 
power and N 0 is the one-sided power spectral density of 
the receiver bandpass noise. 

At lower range of (P/N 0 ) ratio, the receiver may lose 
frequency lock in that the frequency errors at times may 
exceed ±one-half of the sampling frequency, feor ±250 
Hz. Since the error signals e(k) of Equation (11) are the 
same for frequency errors of A Hz as for the case of 
A+nfe Hz for any signed integer n, the estimator may 
make frequency estimation errors of nfe Hz. This situa- 
tion may be referred to as cycle slipping in the fre- 
quency estimator and is akin to the phenomenon of 
cycle slipping (phase errors equal to multiples of 2 77) in 
the phase-estimators. If there are one or more cycles 
slipped in frequency, the computed value of fe, rm5 would 
be much larger compared to the case when no such 
cycle slips occur and would be unacceptable. Thus, 
another important parameter for the performance is the 
probability of maintaining frequency lock throughout 
the trajectory denoted P(lock) or the probability of 
losing the lock Pl= 1 — P(lock). For the purposes of 
estimating the probability by digital computer simula- 
tions, 100 simulation runs are made for each value of 
P/No of interest and an estimate of Pl is plotted vs the 
carrier power to noise power spectral density ratio 
(CNR). The sample rms values of Equation (26) are also 
averaged over all those simulation runs for which the 
frequency lock is maintained. It may well be that for 
sequences that remain under frequency lock there may 
be slipping of cycles in the phase estimates. However, 
even under the presence of such cycle slips, the compu- 
tation made on the basis of Equation (25) provides a 
good estimate of the pseudo-range as evidenced by the 
simulations. One cycle slip only corresponds to an error 
of about 0.2 meters in the pseudo-range estimate. 

FIGS. 6-15 present the simulation results for the 
DLS algorithm and the composite DLS-EKF algo- 
rithms presented in the previous sections. The results 
for the EKF algorithm operating by itself are available 
in V. Vilnrotter, et al., Sept. 1998, supra, for compari- 
son. For the simulation results a value of X equal to 0.97 
has been selected. The initial covariance matrix P(0) for 
the DLS algorithm is selected to be a diagonal matrix 
with its diagonal elements equal to 2X 10 3 , 2X 10 7 and 
2 XlO 9 respectively, reflecting the possible uncertainity 
about the parameters. Three different VCO update 
intervals equal to 5, 10 and 20 sample times have been 
considered. The matrix Q of Equation (19) is also se- 
lected to be a diagonal matrix for convenience, with its 
consecutive diagonal elements equal to 4X 10 3 , 2X 10 6 
and 10 s . The Q matrix represents possible variations in 
the input signal' parameters over an update interval and 
is arrived at from the consideration of a-priori estimate 
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of the maximum possible value of the highest order 
derivative (jerk) present in the input trajectory. For the 
second stage EKF algorithm, the initial covariance 
matrix 2/(0) is selected to be also a diagonal matrix but 
with its elements smaller in value than the correspond- 5 
ing elements of the Q matrix, thus effectively resulting 
in a higher averaging period and smaller estimation 
errors compared to the DLS algorithm. The selected 
values of diagonal elements of 2/0) matrix are equal to 
1.0, 10 3 and 10 6 respectively, in the following Simula- 10 
tions. From the simulations it appeared to be advanta- 
geous in terms of numerical stability to periodically 
reset the covariance matrix P of the DLS algorithm to 
its initial value. Such a period was selected to be 10 
times the VCO update interval. 15 

FIGS. 6 and 7 present the simulation results for the 
performance of DLS algorithm while tracking the high 
dynamic trajectory of FIG. 5 in the absence of any data 
modulation. FIG. 6 plots the probability of losing the 
frequency lock P i as a function of CNR for two differ- 
ent values of N equal to 10 and 20. As may be inferred 
from the figure, a value of Pi of less than 0. 1 is obtained 
for CNRs above 23.1 dB which is defined to be the 
threshold point of the algorithm. FIG. 6 also plots the ^ 
average number of cycles slipped in the frequency esti- 
mation, denoted by N ci . In defining such a cycle slip, 
the entire frequency range is divided into disjoint seg- 
ments of f, Hz with the first segment extending from 
— fs/2 to fj/2 Hz. Whenever the frequency estimation 3Q 
error jumps from one such segment to an adjacent one 
in either direction, a cycle slip is said to occur. FIG. 7 
plots the rms error in the Doppler estimation as com- 
puted from Equation (26) and averaged over all conver- 
gent sequences. For a CNR between 25 and 30 dB-Hz, 35 
an rms error of 10-20 Hz is obtained that corresponds to 
a velocity tracking error of 2-4 m/s. FIG. 7 also plots 
the averaged length L CJ of a slipped cycle Jn terms of 
number of samples. The information about N„and L CT is 
relevant in the case of multistage algorithm. FIGS. 8 
and 9 present the results for the probability of losing 
lock P l and the rms estimation error for the DLS algo- 
rithm in the presence of data modulation for three dif- 
ferent values of N equal to 5, 10 and 20. As may be 
observed from the figures, the presence of data modula- 45 
tion increases the threshold by only 0.25-0.5 dB com- 
pared to the case of no data modulation. The increase in 
rms frequency estimation error is about 10% due to data 
modulation. 

FIGS. 10 and 11 present the performance of the com- 50 
posite DLS-EKF algorithm in the absence of data mod- 
ulation. Note that corresponding to N=5, the threshold 
of the algorithm is 22.75 dB-Hz, which is slightly 
smaller than for the DLS algorithm. However, the rms 
estimation errors are significantly smaller than for the 55 
single-stage DLS algorithm. For the CNR range of 
25-30 dB-Hz, the rms error in the Doppler estimation 
lies in the range of 4-15 Hz corresponding to the veloc- 
ity estimation error range of 0.8 to 3 m/s. The DLS- 
EKF algorithm also provides the carrier phase estimate. 60 
The modulo-27r phase-estimation error is plotted in 
FIG. 12 from which it is clear that the algorithm is 
capable of coherent data detection with small probabil- 
ity of error if the CNR is higher than 25 dB-Hz. In fact, 
as shown in FIG. 15 , the algorithm provides good esti- 65 
mates of pseudo-range (related to the absolute phase 
error via Equation (25)) up to a CNR of about 23 dB at 
which point the rms error is less than 4 m. The rms 
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pseudo range error is less than 1 m for CNRs higher 
than 25.5 dB-Hz. 

The corresponding results for the performance of 
DLS-EKF algorithm in the presence of data modula- 
tion are presented in FIGS. 13 — 15 . For this case, a mini- 
mum threshold of 23.8 dB is obtained for N= 10 which 
is about 1 dB higher than for the case of no data modula- 
tion. In terms of rms estimation errors, for a CNR range 
of about 25-30 dB-Hz, the rms frequency estimation 
error lies in a range of 8-20 Hz corresponding to a 
velocity error of about 1.5 to 4 m/s. For this case, as is 
apparent from FIG. 15 , the pseudo-range estimation 
errors are also higher and for a CNR range of 25-30 
dB-Hz lie in a range of 0.3-6 m. Notice, however, that 
no sharp threshold is observed in either the frequency 
or phase estimation errors over the entire range of CNR 
between 22-30 dB-Hz considered in the simulations. 

7. Comparison with Previous Techniques 

For the case of no data modulation, the following 
compares the performance of the proposed algorithm 
with some of the techniques analyzed in V. Vilnrotter, 
et al., (1988) cited above, in terms of their performance 
when tracking exactly the same high dynamic trajec- 
tory. Compared to a more computation-intensive maxi- 
mum likelihood estimate, the DLS-EKF algorithm re- 
quires about 0.25 dB smaller CNR than MLE in terms 
of threshold. In terms of rms frequency estimation er- 
rors, the MLE achieves an rms error between 8 Hz to 
35 Hz at a CNR of 23 dB-Hz depending upon the esti- 
mation delay ranging between 30-80 samples (higher 
delay provides smaller error). The DLS-EKF algo- 
rithm provides an error of 35 Hz for a delay of 5 sam- 
ples at a CNR of 23 dB-Hz. The MLE algorithm does 
not provide any phase estimate. 

Compared to a single-stage EKF algorithm, the 
DLS-EKF algorithm is bout 1. 5-2.0 dB better in terms 
of threshold depending upon the value of exponential 
data weighting coefficient and the filter order used in 
the simulation of V. Vilnrotter, et al., (1988) cited 
above. In terms of rms errors, the performance is similar 
to that of third-order EKF alone. Notice, however, that 
direct comparison with the results of V. Vilnrotter, et 
al., may be somewhat misleading. This is so because 
while the DLS-EKF algorithm includes all of the se- 
quences in the computation of rms error above 25.5 
dB-Hz, EKF rejects about 5% of the worst sequences, 
as the probability of locking lock is about 0.05 at CNR 
of 25.5 dB-Hz. 

For the case of cross product AFC loop analyzed by 
V. Vilnrotter, et al., the threshold lies in a range of 
25-28 dB-Hz depending upon the loop parameters. 
Thus, the DLS-EKF algorithm is superior by 2-5 
dB-Hz compared to AFC loop. AFC loop provides a 
minimum rms frequency error of 25 Hz at a CNR of 28 
dB-Hz compared to a minimum of 5 Hz achieved for the 
DLS-EKF algorithm for the same CNR. Notice that in 
an AFC loop, the parameters achieving a relatively low 
estimation error are different than those yielding low 
thresholds and thus, a range of loop parameters must be 
considered for proper comparison. In terms of rms 
phase error, the performance of the DLS-EKF algo- 
rithm is similar to EKF alone. In terms of computations, 
the DLS-EKF algorithm requires about the same num- 
ber of computations as for a fourth-order EKF by about 
twice as many computations as a third-order EKF. The 
number of computations are at least an order of magni- 
tude smaller than the MLE. 
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For the case when the data modulation is present, the 
following compares the performance of the DLS-EKF 
algorithm with the MLE algorithm of W. J. Hurd, et al., 
(1987) cited above, where a somewhat less severe GPS 
trajectory is analyzed. The results of W. J. Hurd, et al., 5 
show a marked threshold of about 30 dB-Hz in terms of 
CNR compared to less than 24 dB-Hz for the proposed 
algorithm. Thus, the proposed algorithm results in more 
than 6-dB reduction in terms of threshold compared to 
previous schemes of the literature. 

. In terms of the rms frequency estimation errors at a 
30 dB-Hz CNR, the MLE algorithm provides an rms 
range rate error of 6 m/s compared to an error of less 
than 2 m/s achieved by DLS-EKF algorithm. There is 
also very significant improvement in terms of the rms 15 
position estimation error. At 30 dB-Hz an rms error of 
1 meter is reported by W. S. Hurd, et al., compared to 
about 0.25 meter obtained by the DLS-EKF algorithm. 

In terms of computations, both of the algorithms are 
comparable. In terms of threshold on CNR, the DLS 20 
algorithm is very close to the composite DLS-EKF 
algorithm. However, in terms of rms frequency estima- 
tion errors, it has significantly higher estimation errors. 

In those cases where higher estimation errors are ac- 
ceptable, one may apply the DLS algorithm by itself, as 25 
it requires only one-half of the computations required 
by the DLS-EKF algorithm. 

8. Conclusions 

A novel multistage estimation scheme has been pres- 
ented for the efficient estimation of the phase and fre- 30 
quency of a very high dynamic signal, which may possi- 
bly by phase modulated by unknown binary data and is 
received under relatively low carrier-to-noise power 
ratio conditions. The proposed scheme is of very gen- 
eral nature and has much wider scope than the applica- 35 
tion described as examples. For a very important appli- 
cation of dynamic GPS trajectories, the specific exam- 
ple of a DLS-EKF scheme described has been special- 
ized to have just two stages. The first stage of the esti- 
mation scheme is a least-squares algorithm operating 40 
upon the differential signal model while the second 
stage is an extended Kalman filter of third order. 

For very high dynamic GPS trajectories, the pro- 
posed algorithm has been shown to significantly outper- 
form the previous algorithms reported in the literature 45 
in one or more aspects, including threshold on CNR, 
estimation errors, availability of phase estimates and 
thus, the estimate of pseudo-range, computational com- 
plexity and flexibility. For the case of no data modula- 
tion, the proposed scheme has a threshold that is 50 
slightly lower than the more computation-intensive 
implementation of MLE algorithm. When compared to 
just the EKF operating by itself, the proposed DLS- 
EKF scheme provides from about 1.5 to 2-dB reduction 
in threshold. In comparison to more conventional 55 
schemes, such as AFC loops, the performance is even 
better. 

For the case when an unknown data modulation is 
present, the algorithm provides an improvement of 6 dB 
in terms of threshold on CNR in comparison to the 60 
MLE scheme of W. J. Hurd, et al., (1987) specifically 
proposed for such applications. In addition to phase and 
frequency estimates, the algorithm can provide esti- 
mates of frequency derivative as well, although not 
presented here. The scheme being of a very general 65 
nature, it may be possible to reduce the threshold even 
further by using a higher dimension for the state vector 
related to the higher number of terms in the Taylor 
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series expansion in arriving at the signal model for the 
first stage DLS algorithm'. Further improvements are 
possible by the application of more optimum sampling 
techniques as proposed by R. Kumar, in “Differential 
Sampling for Fast Frequency Acquisition via Adaptive 
Extended Least Squares Algorithm,” (1987) cited 
above. The performance may also be improved both in 
terms of the threshold and the rms estimation errors by 
increasing the number of states to three or more. 

I claim: 

1. A multistage estimator for the parameters of a 
received carrier signal experiencing at least Doppler 
and Doppler rate due to high dynamics in relative mo- 
tion between a receiver and a transmitter of said signal, 
said signal being phase modulated at times by unknown 
data, and where the parameters of interest include 
phase, frequency and frequency derivatives, comprising 
a plurality of estimation stages connected in cascade for 
successively reducing estimation errors due to said dy- 
namics, the first of said stages comprising means for 
estimating parameters of relatively high dynamic sig- 
nals with relatively high estimation errors as compared 
to succeeding stages, each succeeding stage comprising 
means for yielding smaller estimation errors of parame- 
ters of lower dynamic signals, successive ones of said 
estimation stages averaging the received carrier signal 
over progressively longer intervals during which pa- 
rameters of the estimation process remain nearly con- 
stant, thereby progressively reducing the effects of 
noise. 

2. A multistage estimator as defined in claim 1 
wherein said first of said plurality of estimation stages 
connected in cascade is implemented with a differential 
least squares algorithm and the second stage is imple- 
mented with an extended Kalman filter, whereby a 
recursive differential least square and extended Kalman 
filter algorithm is provided for estimation of said param- 
eters of said received signal. 

3. A multistage estimator for the parameters of a 
received carrier signal experiencing at least Doppler 
and Doppler rate due to high dynamics in relative mo- 
tion between a receiver and a transmitter of said signal, 
said signal being phase modulated at times by unknown 
data, where the parameters of interest include phase, 
frequency and frequency derivatives, said estimator 
having a plurality of estimation stages connected in 
cascade for successively reducing estimation errors due 
to said dynamics, a first one of said stages in cascade 
comprising means for using a differential least square 
algorithm for coarse estimation of parameters related to 
frequency of said received signal and its derivatives, but 
with a relatively small probability of the frequency 
estimation error exceeding one-half the sampling fre- 
quency, and a second one of said stages in cascade com- 
prising filter means for operating on the error signal 
available from said first stage to refine the overall esti- 
mates of the phase along with a more refined estimate of 
frequency, and in the process reducing the number of 
cycle slips. 

4. A multistage estimator for parameters of a received 
carrier signal under high dynamic conditions of the 
receiver said estimator having a first stage comprising 
means for coarse estimation of received signal parame- 
ters resulting in higher rms estimation errors but with a 
relatively small probability of the frequency estimation 
error exceeding one half of the sampling frequency, and 
at least a second stage comprising means for operating 
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on the error signals available from said first stage for 
refining the overall estimates of said parameters. 

5. A multistage estimator as defined in claim 4 
wherein said means of said first stage is preferably se- 
lected to use a differential least squares algorithm for 
providing relatively coarse estimates of the received 
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signal frequency and its derivatives, and said means of 
said second stage is preferably selected to be a third- 
order extended Kalman filter which yields a more re- 
fined estimate of frequency as well as an estimate of the 
5 signal phase. 
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